% doing the post processing for the DNS data
clear

iscalar=1;
%Read Kasagi's data
load('kasagi.mat')

yp = importyp('yp.dat');
n_nodes=size(yp,1);

um=read_data('um',n_nodes);
uum=read_data('uum',n_nodes);
uvm=read_data('uvm',n_nodes);
uwm=read_data('uwm',n_nodes);
vm=read_data('vm',n_nodes);
vvm=read_data('vvm',n_nodes);
vwm=read_data('vwm',n_nodes);
wm=read_data('wm',n_nodes);
wwm=read_data('wwm',n_nodes);

Re_m=2280;
Re_tau=150;
u_tau_kasagi=Re_tau/Re_m;
du=dery(yp,um,n_nodes);
du1_m=(abs(du(1))+abs(du(n_nodes)))/2;
if (true)
    u_tau=sqrt(du1_m/Re_m);
else
    u_tau=u_tau_kasagi;
end
h_plus=Re_m*u_tau;
disp(['h_plus ',num2str(h_plus),' is ',num2str((h_plus/Re_tau-1)*100),'% higher than that of Kasagi']);
yp=yp.*h_plus;
um=um./u_tau;
uum=uum./u_tau/u_tau;
uvm=uvm./u_tau/u_tau;
uwm=uwm./u_tau/u_tau;
vm=vm./u_tau;
vvm=vvm./u_tau/u_tau;
vwm=vwm./u_tau/u_tau;
wm=wm./u_tau;
wwm=wwm./u_tau/u_tau;
    
x_range=[0 Re_tau];
linear_log='linear';
%linear_log='log';

Ruu=sqrt(uum-um.*um);
Ruv=sqrt(abs(uvm-um.*vm));
Ruw=sqrt(abs(uwm-um.*wm));
Rvv=sqrt(vvm-vm.*vm);
Rvw=sqrt(abs(vwm-vm.*wm));
Rww=sqrt(wwm-wm.*wm);

if (true)
    figure1=figure;
    plot(kasagi.y,kasagi.u_mean);
    hold('all');
    plot(yp,um,'-.');
% function figure_property( figure_control, figure_name, x_label, y_label, x_range, linear_log)
% e.g.  x_range=[0 110]
% linear_log='log'
    figure_property( figure1, 'u\_mean vs y+', 'y+', 'u\_mean', x_range, linear_log );

    figure2=figure;
    plot(kasagi.y,kasagi.uu);
    hold('all');
    plot(yp,Ruu,'-.');
    figure_property( figure2, 'u\_rms vs y+', 'y+', 'u\_rms', x_range, linear_log );

    figure3=figure;
    plot(kasagi.y,kasagi.vv);
    hold('all');
    plot(yp,Rvv,'-.');
    figure_property( figure3, 'v\_rms vs y+', 'y+', 'v\_rms', x_range, linear_log );

    figure4=figure;
    plot(kasagi.y,kasagi.ww);
    hold('all');
    plot(yp,Rww,'-.');
    figure_property( figure4, 'w\_rms vs y+', 'y+', 'w\_rms', x_range, linear_log );
end

if (iscalar==1)
    Pr=0.025;
    tm=read_data('phim',n_nodes);
    ttm=read_data('phiphim',n_nodes);
    utm=read_data('uphim',n_nodes);
    vtm=read_data('vphim',n_nodes);
    wtm=read_data('wphim',n_nodes);
    dt=dery(yp,tm,n_nodes);
    dt1_m=(abs(dt(1))+abs(dt(n_nodes)))/2;
    t_tau=sqrt(dt1_m/h_plus/Pr);
    tm=tm./t_tau;
    ttm=ttm./t_tau/t_tau;
    utm=utm./t_tau/u_tau;
    vtm=vtm./t_tau/u_tau;
    wtm=wtm./t_tau/u_tau;
    t_rms=sqrt(ttm-tm.*tm);
    
    
    figure5=figure;
    plot(kasagi.y,kasagi.t_mean);
    hold('all');
    plot(yp,tm,'-.');
    figure_property( figure5, 't\_mean vs y+', 'y+', 't_mean', x_range, linear_log );

    figure6=figure;
    plot(kasagi.y,kasagi.t_rms);
    hold('all');
    plot(yp,t_rms,'-.');
    figure_property( figure6, 't\_rms vs y+', 'y+', 't\_rms', x_range, linear_log );

    figure7=figure;
    plot(kasagi.y,kasagi.vt);
    hold('all');
    plot(yp,-vtm,'-.');
    figure_property( figure7, '< vt > vs y+', 'y+', '< vt >', x_range, linear_log );

end